home *** CD-ROM | disk | FTP | other *** search
/ CU Amiga Super CD-ROM 22 / CU Amiga Magazine's Super CD-ROM 22 (1998)(EMAP Images)(GB)[!][issue 1998-05].iso / PowerPC / Programming / PPCsiod / sources / math.c < prev    next >
Encoding:
C/C++ Source or Header  |  1993-09-25  |  8.8 KB  |  356 lines
  1. /* Scheme In One Define.
  2.  
  3. The garbage collector, the name and other parts of this program are
  4.  
  5.  *                     COPYRIGHT (c) 1989 BY                              *
  6.  *      PARADIGM ASSOCIATES INCORPORATED, CAMBRIDGE, MASSACHUSETTS.       *
  7.  
  8. Conversion  to  full scheme standard, characters, vectors, ports, complex &
  9. rational numbers, and other major enhancments by
  10.  
  11.  *      Scaglione Ermanno, v. Pirinoli 16 IMPERIA P.M. 18100 ITALY        * 
  12.  
  13. Permission  to use, copy, modify, distribute and sell this software and its
  14. documentation  for  any purpose and without fee is hereby granted, provided
  15. that  the  above  copyright  notice appear in all copies and that both that
  16. copyright   notice   and   this  permission  notice  appear  in  supporting
  17. documentation,  and that the name of Paradigm Associates Inc not be used in
  18. advertising or publicity pertaining to distribution of the software without
  19. specific, written prior permission.
  20.  
  21. PARADIGM  DISCLAIMS  ALL WARRANTIES WITH REGARD TO THIS SOFTWARE, INCLUDING
  22. ALL  IMPLIED  WARRANTIES  OF MERCHANTABILITY AND FITNESS, IN NO EVENT SHALL
  23. PARADIGM  BE  LIABLE  FOR ANY SPECIAL, INDIRECT OR CONSEQUENTIAL DAMAGES OR
  24. ANY DAMAGES WHATSOEVER RESULTING FROM LOSS OF USE, DATA OR PROFITS, WHETHER
  25. IN  AN ACTION OF CONTRACT, NEGLIGENCE OR OTHER TORTIOUS ACTION, ARISING OUT
  26. OF OR IN CONNECTION WITH THE USE OR PERFORMANCE OF THIS SOFTWARE.
  27.  
  28. */
  29.  
  30. #include <stdio.h>
  31. #include <string.h>
  32. #include <ctype.h>
  33. #include <setjmp.h>
  34. #include <signal.h>
  35. #include <math.h>
  36. #include <limits.h>
  37.  
  38. #include "siod.h"
  39.  
  40. #define fixhash(f)((((f)+2357)*1193)%fixarray_dim)
  41.  
  42. LISP intcons(long x)
  43. {LISP z,l;
  44.  long hash,flag;
  45.  hash = fixhash(x);
  46.  for(l=fixarray[hash];CONSP(l);l=CDR(l))
  47.    if (INTNM(CAR(l)) == x)
  48.       return(CAR(l));
  49.  flag=no_interrupt(1);
  50.  NEWCELL(z,tc_intnum);
  51.  INTNM(z)=x;
  52.  fixarray[hash] = cons(z,fixarray[hash]);
  53.  no_interrupt(0);
  54.  return(z);}
  55.  
  56. LISP ratcons(double x,double y)
  57. {double gcd,integ;
  58.  if(y == 0.) err("division by zero",NIL,ERR_GEN);
  59.  if(x == 0.) return(intcons(0));
  60.  if((modf(x,&integ) != 0.) || (modf(y,&integ) != 0.))
  61.    return flocons(x/y);
  62.  gcd = cal_gcd_double(x,y);
  63.  x /= gcd;
  64.  y /= gcd;
  65.  if(y<0.)
  66.   {y=-y;
  67.    x=-x;}
  68.  if(y == 1.)
  69.    return flocons(x);
  70.  if((x>LONG_MAX) || (x<LONG_MIN) || (y>ULONG_MAX))
  71.    return flocons(x/y);
  72.  return(ratco((long)x,(unsigned long)y));}
  73.  
  74. LISP ratco(long x,unsigned long y)
  75. {long flag;
  76.  LISP z;
  77.  flag=no_interrupt(1);
  78.  NEWCELL(z,tc_ratnum);
  79.  RATNUM(z) = x;
  80.  RATDEN(z) = y;
  81.  no_interrupt(flag);
  82.   return(z);}
  83.  
  84. LISP flocons(double x)
  85. {double integ,frac;
  86.  frac = modf(x,&integ); 
  87.  if((frac==0.) && (integ<LONG_MAX) && (integ>LONG_MIN))
  88.    return(intcons((long)x));
  89.  return(floco(x));}
  90.  
  91. LISP floco(double x)
  92. {long flag;
  93.  LISP z;
  94.  flag=no_interrupt(1);
  95.  NEWCELL(z,tc_flonum);
  96.  FLONM(z) = x;
  97.  no_interrupt(flag);
  98.   return(z);}
  99.  
  100. LISP compcons(float x,float y)
  101. {if(y == 0.) return(flocons((double)x));
  102.  return(compco(x,y));}
  103.  
  104. LISP compco(float x,float y)
  105. {long flag;
  106.  LISP z;
  107.  flag=no_interrupt(1);
  108.  NEWCELL(z,tc_compnum);
  109.  COMPRE(z) = x;
  110.  COMPIM(z) = y;
  111.  no_interrupt(flag);
  112.   return(z);}
  113.  
  114. LISP tofloat(LISP x)
  115. {LISP z;
  116.  switch(TYPE(x))
  117.  {case tc_intnum:
  118.     z = floco((double)INTNM(x));
  119.     break;
  120.   case tc_ratnum:
  121.     z = floco((double)RATNUM(x)/(double)RATDEN(x));
  122.     break;
  123.   case tc_flonum:
  124.     z = floco(FLONM(x));
  125.     break;   
  126.   case tc_compnum:
  127.     if(COMPIM(x)==0.)
  128.       z = floco((double)COMPRE(x));
  129.     else
  130.       z = compcons(COMPRE(x),COMPIM(x));
  131.     break;}
  132.  return(z);}
  133.  
  134. LISP ltofloat(LISP x)
  135. {if(NNUMBERP(x))
  136.   err("float",x,ERR_GEN_ARG | ERR_NNUM);
  137.  return(tofloat(x));}
  138.  
  139. LISP torational(LISP x)
  140. {LISP z;
  141.  double integ;
  142.  switch(TYPE(x))
  143.  {case tc_intnum:
  144.     z = x;
  145.     break;
  146.   case tc_ratnum:
  147.     z = ratco(RATNUM(x),RATDEN(x)); 
  148.     break;
  149.   case tc_flonum:
  150.     modf((FLONM(x) * 100000.),&integ);
  151.     z = ratcons(integ,100000.);
  152.     break;   
  153.   case tc_compnum:
  154.     if(COMPIM(x)==0.)
  155.      {modf((COMPRE(x) * 100000.),&integ);
  156.       z = ratcons(integ,100000.);}
  157.     else
  158.       z = compcons(COMPRE(x),COMPIM(x));
  159.     break;}
  160.  if(INTNUMP(z))
  161.     z = ratco(INTNM(x),1);
  162.  return(z);}
  163.  
  164. LISP ltorational(LISP x)
  165. {if(NNUMBERP(x))
  166.   err("rational",x,ERR_GEN_ARG | ERR_NNUM);
  167.  return(torational(x));}
  168.  
  169. LISP tocomplex(LISP x)
  170. {LISP z;
  171.  switch(TYPE(x))
  172.  {case tc_intnum:
  173.     z = compco((float)INTNM(x),0.F); 
  174.     break;
  175.   case tc_ratnum:
  176.     z = compco((float)((float)RATNUM(x)/(float)RATDEN(x)),0.F);
  177.     break;
  178.   case tc_flonum:
  179.     z = compco((float)FLONM(x),0.F);
  180.     break;
  181.   case tc_compnum:
  182.     z = compco(COMPRE(x),COMPIM(x));
  183.     break;}
  184.  return(z);}
  185.  
  186. LISP ltocomplex(LISP x)
  187. {if(NNUMBERP(x))
  188.   err("complex",x,ERR_GEN_ARG | ERR_NNUM);
  189.  return(tocomplex(x));}
  190.  
  191. LISP converti(LISP x,short y)
  192. {LISP z;
  193.  switch(y)
  194.   {case tc_intnum:
  195.      z = ltruncate(x);
  196.      break;
  197.    case tc_ratnum:
  198.      z = torational(x);
  199.      break;
  200.    case tc_flonum:
  201.      z = tofloat(x);
  202.      break;
  203.    case tc_compnum:
  204.      z = tocomplex(x);
  205.      break;} 
  206.  return(z);
  207. }
  208.  
  209. LISP plus2(LISP x,LISP y)
  210. {LISP z;
  211.  if(TYPE(x) > TYPE(y))
  212.      y=converti(y,TYPE(x));
  213.  else if(TYPE(x) < TYPE(y))
  214.      x=converti(x,TYPE(y));
  215.  switch(TYPE(x))
  216.  {
  217.   case tc_intnum:
  218.     z = flocons(((double)INTNM(x)+(double)INTNM(y)));
  219.     break;
  220.   case tc_ratnum:
  221.     z = ratcons((((double)RATNUM(x) * 
  222.                         (double)RATDEN(y)) +
  223.                        ((double)RATNUM(y) * 
  224.                         (double)RATDEN(x))),
  225.                 ((double)RATDEN(x)*(double)RATDEN(y)));
  226.     break;
  227.   case tc_flonum:
  228.     z = flocons(FLONM(x)+FLONM(y));
  229.     break;
  230.   case tc_compnum:
  231.     z = compcons((float)(COMPRE(x)+COMPRE(y)),(float)(COMPIM(x)+COMPIM(y)));
  232.     break;}
  233.   return(z);}
  234.  
  235. LISP minus2(LISP x,LISP y)
  236. {LISP z;
  237.  if(TYPE(x) > TYPE(y))
  238.      y=converti(y,TYPE(x));
  239.  else if(TYPE(x) < TYPE(y))
  240.      x=converti(x,TYPE(y));
  241.  switch(TYPE(x))
  242.  {
  243.   case tc_intnum:
  244.     z = flocons(((double)INTNM(x)-(double)INTNM(y)));
  245.     break;
  246.   case tc_ratnum:
  247.     z = ratcons((((double)RATNUM(x) * 
  248.                         (double)RATDEN(y)) -
  249.                        ((double)RATNUM(y) * 
  250.                         (double)RATDEN(x))),
  251.                 ((double)RATDEN(x)*(double)RATDEN(y)));
  252.     break;
  253.   case tc_flonum:
  254.     z = flocons(FLONM(x)-FLONM(y));
  255.     break;
  256.   case tc_compnum:
  257.     z = compcons((float)(COMPRE(x)-COMPRE(y)),(float)(COMPIM(x)-COMPIM(y)));
  258.     break;}
  259.   return(z);}
  260.  
  261. LISP times2(LISP x,LISP y)
  262. {LISP z;
  263.  if(TYPE(x) > TYPE(y))
  264.      y=converti(y,TYPE(x));
  265.  else if(TYPE(x) < TYPE(y))
  266.      x=converti(x,TYPE(y));
  267.  switch(TYPE(x))
  268.  {
  269.   case tc_intnum:
  270.     z = flocons(((double)INTNM(x)*(double)INTNM(y)));
  271.     break;
  272.   case tc_ratnum:
  273.     z = ratcons(((double)RATNUM(x)*(double)RATNUM(y)),
  274.                 ((double)RATDEN(x)*(double)RATDEN(y)));
  275.     break;
  276.   case tc_flonum:
  277.     z = flocons(FLONM(x)*FLONM(y));
  278.     break;
  279.   case tc_compnum:
  280.     z = compcons((float)((COMPRE(x)*COMPRE(y))-(COMPIM(x)*COMPIM(y))),
  281.                  (float)((COMPRE(x)*COMPIM(y))+(COMPRE(y)*COMPIM(x))));
  282.     break;}
  283.   return(z);}
  284.  
  285. LISP divide2(LISP x,LISP y)
  286. {LISP z;
  287.  if(TYPE(x) > TYPE(y))
  288.      y=converti(y,TYPE(x));
  289.  else if(TYPE(x) < TYPE(y))
  290.      x=converti(x,TYPE(y));
  291.  switch(TYPE(x))
  292.  {
  293.   case tc_intnum:
  294.        z = ratcons((double)INTNM(x),(double)INTNM(y));
  295.     break;
  296.   case tc_ratnum:
  297.        z = ratcons(((double)RATNUM(x)*(double)RATDEN(y)),
  298.                    ((double)RATDEN(x)*fabs((double)RATNUM(y))));
  299.     break;
  300.   case tc_flonum:
  301.     z = flocons(FLONM(x)/FLONM(y));
  302.     break;
  303.   case tc_compnum:
  304.    {float mo = (COMPRE(y)*COMPRE(y))+(COMPIM(y)*COMPIM(y));
  305.     z = compcons((float)(((COMPRE(x)*COMPRE(y))+(COMPIM(x)*COMPIM(y)))/mo),
  306.                  (float)(((COMPRE(y)*COMPIM(x))-(COMPRE(x)*COMPIM(y)))/mo));
  307.     break;}}
  308.   return(z);}
  309.  
  310. LISP getnumer(LISP num)
  311. {LISP z;
  312.  if(NRATNUMP(num))err("arg to numerator must be a rational",num,ERR_GEN);
  313.  z = intcons(RATNUM(num));
  314.  return(z);}
  315.  
  316. LISP getdenom(LISP num)
  317. {LISP z;
  318.  if(NRATNUMP(num))err("arg to denominator must be a rational",num,ERR_GEN);
  319.  z = intcons(RATDEN(num));
  320.  return(z);}
  321.  
  322. LISP makerat(LISP num,LISP den)
  323. {LISP z;
  324.  double dum;
  325.  num = tofloat(num);
  326.  den = tofloat(den);
  327.  if(NFLONUMP(num)||(modf(FLONM(num),&dum)!=0.))
  328.     err("make-rational",num,ERR_GEN_ARG | ERR_NINT);
  329.  if(NFLONUMP(den)||(modf(FLONM(den),&dum)!=0.))
  330.     err("make-rational",num,ERR_GEN_ARG | ERR_NINT);
  331.  z = ratcons(FLONM(num),FLONM(den));
  332.  return(z);}
  333.  
  334. LISP makecomp(LISP real,LISP imag)
  335. {LISP z;
  336.  real = tofloat(real);
  337.  imag = tofloat(imag);
  338.  if(NFLONUMP(real))
  339.     err("arg to make-complex must be a float",real,ERR_GEN);
  340.  if(NFLONUMP(imag)) 
  341.     err("arg to make-complex must be a float",imag,ERR_GEN);
  342.  z = compcons((float)FLONM(real),(float)FLONM(imag));
  343.  return(z);}
  344.  
  345. LISP getreal(LISP num)
  346. {LISP z;
  347.  if(NCOMPNUMP(num))err("arg to real must be a complex",num,ERR_GEN);
  348.  z = flocons((double)COMPRE(num));
  349.  return(z);}
  350.  
  351. LISP getimag(LISP num)
  352. {LISP z;
  353.  if(NCOMPNUMP(num))err("arg to imaginary must be a complex",num,ERR_GEN);
  354.  z = flocons((double)COMPIM(num));
  355.  return(z);}
  356.